function results = correlated_mimo_capacity(varargin)
% 相关信道中的MIMO容量分析
% 输入参数:
%   'SNR_dB' - SNR范围 (dB), 默认: -5:5:30
%   'correlation_values' - 相关性值数组, 默认: [0, 0.3, 0.6, 0.9]
%   'Nt', 'Nr' - 天线数量, 默认: 4x4
%   'num_realizations' - 蒙特卡洛仿真次数, 默认: 500
% 输出:
%   results - 包含容量结果、特征值分布、条件数等信息的结构体

% 解析输入参数
p = inputParser;
addParameter(p, 'SNR_dB', -5:5:30);
addParameter(p, 'correlation_values', [0, 0.3, 0.6, 0.9]);
addParameter(p, 'Nt', 4);
addParameter(p, 'Nr', 4);
addParameter(p, 'num_realizations', 500);
parse(p, varargin{:});

SNR_dB = p.Results.SNR_dB;
correlation_values = p.Results.correlation_values;
Nt = p.Results.Nt;
Nr = p.Results.Nr;
num_realizations = p.Results.num_realizations;
SNR_linear = 10.^(SNR_dB/10);

% 添加路径
addpath('../Common');

% 获取颜色定义
colors = color_definitions();

fprintf('=== 相关信道中的MIMO容量分析 ===\n');

% 空间相关性影响
corr_capacity = zeros(length(SNR_dB), length(correlation_values));

for corr_idx = 1:length(correlation_values)
    rho = correlation_values(corr_idx);
    fprintf('空间相关性: %.1f\n', rho);
    
    for snr_idx = 1:length(SNR_dB)
        snr = SNR_linear(snr_idx);
        
        % 生成相关信道
        capacities_corr = zeros(num_realizations, 1);
        
        for real = 1:num_realizations
            H_corr = generate_correlated_channel(Nt, Nr, rho);
            capacities_corr(real) = log2(det(eye(Nr) + (snr/Nt) * H_corr * H_corr'));
        end
        
        corr_capacity(snr_idx, corr_idx) = mean(capacities_corr);
    end
end

% 绘制相关信道容量
figure('Name', '相关信道中的MIMO容量', 'Position', [250, 250, 1200, 800]);

subplot(2,2,1);
for corr_idx = 1:length(correlation_values)
    plot(SNR_dB, corr_capacity(:, corr_idx), ...
         ['-', colors(corr_idx)], 'LineWidth', 2);
    hold on;
end
grid on;
xlabel('SNR (dB)');
ylabel('容量 (bps/Hz)');
title('相关信道中的MIMO容量');
legend(arrayfun(@(x) sprintf('ρ=%.1f', x), correlation_values, 'UniformOutput', false));

% 容量损失分析
subplot(2,2,2);
capacity_loss = (corr_capacity(:, 1) - corr_capacity) ./ corr_capacity(:, 1) * 100;
for corr_idx = 2:length(correlation_values)
    plot(SNR_dB, capacity_loss(:, corr_idx), ...
         ['-', colors(corr_idx)], 'LineWidth', 2);
    hold on;
end
grid on;
xlabel('SNR (dB)');
ylabel('容量损失 (%)');
title('相关性导致的容量损失');
legend(arrayfun(@(x) sprintf('ρ=%.1f', x), correlation_values(2:end), 'UniformOutput', false));

% 特征值分布分析 (固定SNR)
snr_fixed = 15; % dB
snr_idx = find(SNR_dB == snr_fixed, 1);

eigenvalue_distributions = cell(length(correlation_values), 1);
for corr_idx = 1:length(correlation_values)
    rho = correlation_values(corr_idx);
    
    % 生成多个信道样本
    num_samples = 1000;
    all_eigenvalues = [];
    
    for i = 1:num_samples
        H_sample = generate_correlated_channel(Nt, Nr, rho);
        eigenvals = eig(H_sample * H_sample');
        all_eigenvalues = [all_eigenvalues; eigenvals];
    end
    
    eigenvalue_distributions{corr_idx} = all_eigenvalues;
end

subplot(2,2,3);
for corr_idx = 1:length(correlation_values)
    histogram(eigenvalue_distributions{corr_idx}, 30, 'Normalization', 'pdf');
    hold on;
end
grid on;
xlabel('特征值');
ylabel('概率密度');
title(sprintf('特征值分布 (SNR=%d dB)', snr_fixed));
legend(arrayfun(@(x) sprintf('ρ=%.1f', x), correlation_values, 'UniformOutput', false));

% 条件数统计
subplot(2,2,4);
condition_numbers = zeros(length(correlation_values), 1);
for corr_idx = 1:length(correlation_values)
    rho = correlation_values(corr_idx);
    
    cond_nums = zeros(num_realizations, 1);
    for real = 1:num_realizations
        H_cond = generate_correlated_channel(Nt, Nr, rho);
        cond_nums(real) = cond(H_cond);
    end
    
    condition_numbers(corr_idx) = mean(cond_nums);
end

bar(condition_numbers);
grid on;
title('信道条件数 vs 相关性');
xlabel('相关性');
ylabel('平均条件数');
set(gca, 'XTickLabel', arrayfun(@(x) sprintf('%.1f', x), correlation_values, 'UniformOutput', false));

% 组织结果
results.SNR_dB = SNR_dB;
results.correlation_values = correlation_values;
results.corr_capacity = corr_capacity;
results.capacity_loss = capacity_loss;
results.eigenvalue_distributions = eigenvalue_distributions;
results.condition_numbers = condition_numbers;
results.Nt = Nt;
results.Nr = Nr;

end

function H_corr = generate_correlated_channel(Nt, Nr, rho)
    % 生成相关MIMO信道
    % 发送端相关矩阵
    R_tx = ones(Nt, Nt);
    for i = 1:Nt
        for j = 1:Nt
            R_tx(i,j) = rho^abs(i-j);
        end
    end
    
    % 接收端相关矩阵
    R_rx = ones(Nr, Nr);
    for i = 1:Nr
        for j = 1:Nr
            R_rx(i,j) = rho^abs(i-j);
        end
    end
    
    % 相关信道
    R_sqrt_tx = sqrtm(R_tx);
    R_sqrt_rx = sqrtm(R_rx);
    H_uncorr = sqrt(0.5) * (randn(Nr, Nt) + 1i * randn(Nr, Nt));
    H_corr = R_sqrt_rx * H_uncorr * R_sqrt_tx';
end